/*!
 * \file OsterCliquer.cpp
 *
 * \date Jan 27, 2014
 * \author correa
 */

#include <cstring>
#include <fstream>
#include <cstdlib>
#include <OsterCliquer.h>

OsterCliquer::OsterCliquer(const Grafo &g, bool clqOrStab) :
		_n(g.getN()), _grafo(g), _RR(clqOrStab ? *(new VertexSet(_n)) : getComplParamSet(_n)),
		_WW(_n), _nodeR(), _nodeRR(), _nodeW(),
		_stack(new Doll*[_n]), _C(new double[_n])
{
	if (clqOrStab)
	{
		for (int i = 0; i < _n; i++)
			_stack[i] = new Doll(_grafo, new VertexSet(_n));
	}
	else
	{
		for (int i = 0; i < _n; i++)
			_stack[i] = new Doll(_grafo, &getComplParamSet(_n));
	}
}

OsterCliquer::~OsterCliquer()
{
	delete &_RR;

	for (int i = 0; i < _n; i++)
		delete _stack[i];

	delete[] _stack;
	delete[] _C;
}

/*!
 * \fn int greedyColoring(VertexSet& R, VertexSet& W, int lim)
 * \brief Greedy heuristic for vertex coloring
 *
 * The vertices in the specified bitmap are colored, removed from this bitmap, and stored in another specified bitmap.
 * A maximum number of colors is specified. Vertices that cannot be colored with this number of colors remains in the input bitmap.
 *
 * \param R Input bitmap of vertices to color.
 * \param W Output bitmap to store vertices as soon as they are colored.
 * \param lim Maximum number of colors that can be used.
 *
 * \return The number of colors used.
 */
int OsterCliquer::greedyColoring(VertexSet& R, VertexSet& W, int lim) {
	int size = 0;
	int vidx;
	int nodeidx = 0;
	bool nemptcol = true;
	int lastnode = NODEIDX(_n-1);

	while (size < lim && nemptcol) {
		_RR.copyAll(R);
		nodeidx = 0;
		nemptcol = false;
		while (nodeidx <= lastnode) {
			_RR.nodeView(nodeidx, _nodeRR);
			R.nodeView(nodeidx, _nodeR);
			W.nodeView(nodeidx, _nodeW);
			vidx = _nodeRR.least();
			while (vidx >= 0) {
				nemptcol = true;
				// note that vidx is removed from _RR only if clqOrStab=false
				_RR.removeAll(_grafo.vecinos(NELEM(nodeidx)+vidx));
				_nodeRR.remove(vidx);
				_nodeR.remove(vidx);
				_nodeW.add(vidx);
				vidx = _nodeRR.least();
			}
			nodeidx++;
		}

		if (nemptcol)
			size++;
	}

	return size;
}

int OsterCliquer::greedyColoring(VertexSet& R) {
	int size = 0;
	int vidx;
	int nodeidx = 0;
	bool nemptcol = true;
	int lastnode = NODEIDX(_n-1);

	_WW.removeAll();

	while (nemptcol) {
		_RR.copyAll(R);
		nodeidx = 0;
		nemptcol = false;
		while (nodeidx <= lastnode) {
			_RR.nodeView(nodeidx, _nodeRR);
			R.nodeView(nodeidx, _nodeR);
			_WW.nodeView(nodeidx, _nodeW);
			vidx = _nodeRR.least();
			while (vidx >= 0) {
				nemptcol = true;
				// note that vidx is removed from _RR only if clqOrStab=false
				_RR.removeAll(_grafo.vecinos(NELEM(nodeidx)+vidx));
				_nodeRR.remove(vidx);
				_nodeR.remove(vidx);
				_nodeW.add(vidx);
				vidx = _nodeRR.least();
			}
			nodeidx++;
		}

		if (nemptcol)
			size++;
	}

	R.addAll(_WW);

	return size;
}

/*!
 * \fn void maxStab(VertexSet& R, VertexSet& W, int u)
 * \brief Greedy heuristic for maximum stable set containing a specified vertex.
 *
 * The vertices in the specified bitmap are greedily searched for a maximum stable set. The vertices in the stable set found are removed from
 * the input bitmap and stored in the output bitmap.
 *
 * \param R Input bitmap of vertices to color.
 * \param W Output bitmap to store vertices as soon as they are colored.
 * \param u The vertex to expand to a stable set.
 */
void OsterCliquer::maxStab(VertexSet& R, VertexSet& W, int u) {
	int nodeidx = 0;
	_RR.copyAll(R);

	// note that u is removed from _RR only if clqOrStab=false
	// usually, u is not in _RR
	_RR.removeAll(_grafo.vecinos(u));

	int lastnode = NODEIDX(_n-1);
	while (_RR.nodeView(nodeidx, _nodeRR).isEmpty() && ++nodeidx <= lastnode);

	while (nodeidx <= lastnode) {
		_RR.nodeView(nodeidx, _nodeRR);
		R.nodeView(nodeidx, _nodeR);
		W.nodeView(nodeidx, _nodeW);
		int vidx = _nodeRR.least();
		while (vidx >= 0) {
			// note that vidx is removed from _RR only if clqOrStab=false
			_RR.removeAll(_grafo.vecinos(NELEM(nodeidx)+vidx));
			_nodeRR.remove(vidx);
			_nodeR.remove(vidx);
			_nodeW.add(vidx);
			vidx = _nodeRR.least();
		}
		nodeidx++;
	}
}

/*!
 * \fn int expClique(VertexSet& R, int& C[], int topC)
 * \brief Greedy heuristic for expanding a specified clique.
 *
 * The vertices in the specified bitmap are greedily searched to entering the specified clique.
 *
 * \param R Input bitmap of vertices to color.
 * \param C The input clique.
 * \param topC Size of input clique.
 *
 * \return The size of the expanded clique.
 */
int OsterCliquer::expClique(VertexSet& R, int * C, int * topC) {
	_RR.copyAll(R);

	int size = 0;

	int t = *topC;
	while (--t >= 0) {
		// note that vidx only remains in _RR if clqOrStab=false
		_RR.retainAll(_grafo.vecinos(C[t]));
		size++;
	}

	if (_RR.isEmpty())
		return *topC;

	int lastnode = NODEIDX(_n-1);
	int nodeidx = 0;
	int v;
	while (nodeidx <= lastnode) {
		_RR.nodeView(nodeidx, _nodeRR);
		int vidx = _nodeRR.least();
		while (vidx >= 0) {
			v = NELEM(nodeidx)+vidx;
			C[size++] = v;
			// note that vidx only remains in _RR if clqOrStab=false
			_RR.retainAll(_grafo.vecinos(v));
			_nodeRR.remove(vidx);
			vidx = _nodeRR.least();
		}
		nodeidx++;
	}

	return size;
}

int OsterCliquer::optimize(const VertexSet& subgrafo, VertexSet& clique, int cotaInferior, int cotaSuperior,
		int tolerancia)
{
	_stack[0]->_R->copyAll(subgrafo);
	_stack[0]->nodeidx = 0;
	_stack[0]->_R->nodeView(_stack[0]->nodeidx, *_stack[0]->Rnode);
	_stack[0]->_W->removeAll();

	int lastnode = NODEIDX(_n-1);

	int curC[_n];		//!< Current clique
	int topC;
	Doll * doll;
	int	v;
	int submax;

	//initialize with a maximal clique
	int expsz = 0;
	int i = expClique(*_stack[0]->_R, curC, &expsz);
	if (i > cotaInferior) {
		cotaInferior = i;
		clique.removeAll();
		for (int h = 0; h < expsz; h++)
			clique.add(curC[h]);

		int ub = greedyColoring(*_stack[0]->_R);
		if (cotaSuperior > cotaInferior + ub)
			cotaSuperior = cotaInferior + ub;

		if (cotaSuperior-cotaInferior <= tolerancia)
			goto TERMINATE;
	}

	topC = 0;
	doll = _stack[topC];
	v = doll->nextUpTo(lastnode);
	submax = 0;
	while (v >= 0 && submax < cotaInferior)
	{
		submax++;
		_C[v] = submax;

		v = doll->nextUpTo(lastnode);
	}

	if (v < 0)
		goto TERMINATE;

	curC[topC++] = v;
	_stack[topC]->startUp(v, *_stack[topC-1]->_W);
	doll = _stack[topC];

	while (topC >= 0 && cotaSuperior-cotaInferior > tolerancia) {
		if (topC == 0)
			_C[curC[topC]] = submax;
		if (topC > 0 && doll->_R->isEmpty())
		{
			if (topC > submax)
			{
				submax = topC;

				if (submax > cotaInferior) {
					submax = expClique(*_stack[0]->_R, curC, &submax);
					clique.removeAll();
					for (int h = 0; h < submax; h++)
						clique.add(curC[h]);

					_stack[0]->_R->add(curC[0]);
					int ub = greedyColoring(*_stack[0]->_R);
					if (cotaSuperior > cotaInferior + ub)
						cotaSuperior = cotaInferior + ub;
					_stack[0]->_R->remove(curC[0]);

					v = _stack[0]->nextUpTo(lastnode);
					_C[v] = (_grafo.isArista(v, curC[0]) ? ++cotaInferior : cotaInferior);
					while (cotaInferior < submax)
					{
						v = _stack[0]->nextUpTo(lastnode);
						_C[v] = ++cotaInferior;
					}
				}
			}

			topC--;
			if (topC >= 0)
				doll = _stack[topC];
		}
		else
		{
			v = doll->nextUpTo(lastnode);
			if (v < 0) {
				topC--;
				if (topC >= 0)
					doll = _stack[topC];
			}
			else if (topC == 0 || topC + _C[v] > cotaInferior)
			{
				curC[topC++] = v;

				//Calculate the new candidate vector
				_stack[topC]->startUp(v, *doll->_W);
				doll = _stack[topC];
			}
			else if (topC + _C[v] > submax)
				submax = topC + _C[v];
		}
	}

	TERMINATE:

	return cotaInferior;
}

int OsterCliquer::greedyColoring(VertexSet& R, const IntValueArray& pesos) {
	int weight = 0, cw;
	int vidx;
	int nodeidx = 0;
	bool nemptcol = true;
	int lastnode = NODEIDX(_n-1);

	_WW.removeAll();

	while (nemptcol) {
		_RR.copyAll(R);
		nodeidx = 0;
		nemptcol = false;
		cw = 0;
		while (nodeidx <= lastnode) {
			_RR.nodeView(nodeidx, _nodeRR);
			R.nodeView(nodeidx, _nodeR);
			_WW.nodeView(nodeidx, _nodeW);
			vidx = _nodeRR.least();
			while (vidx >= 0) {
				int v = NELEM(nodeidx)+vidx;
				nemptcol = true;
				if (pesos[v] > cw)
					cw = pesos[v];
				// note that vidx is removed from _RR only if clqOrStab=false
				_RR.removeAll(_grafo.vecinos(v));
				_nodeRR.remove(vidx);
				_nodeR.remove(vidx);
				_nodeW.add(vidx);
				vidx = _nodeRR.least();
			}
			nodeidx++;
		}

		if (nemptcol)
			weight+=cw;
	}

	R.addAll(_WW);

	return weight;
}

int OsterCliquer::expClique(VertexSet& R, int * C, int * topC, const IntValueArray& pesos) {
	_RR.copyAll(R);

	int weight = 0;

	int t = *topC;
	while (--t >= 0) {
		// note that vidx only remains in _RR if clqOrStab=false
		_RR.retainAll(_grafo.vecinos(C[t]));
		weight+=pesos[C[t]];
	}

	if (_RR.isEmpty())
		return weight;

	int lastnode = NODEIDX(_n-1);
	int nodeidx = 0;
	int v;
	while (nodeidx <= lastnode) {
		_RR.nodeView(nodeidx, _nodeRR);
		int vidx = _nodeRR.least();
		while (vidx >= 0) {
			v = NELEM(nodeidx)+vidx;
			C[(*topC)++] = v;
			weight += pesos[v];
			// note that vidx only remains in _RR if clqOrStab=false
			_RR.retainAll(_grafo.vecinos(v));
			_nodeRR.remove(vidx);
			vidx = _nodeRR.least();
		}
		nodeidx++;
	}

	return weight;
}

int OsterCliquer::optimize(const IntValueArray& pesos, const VertexSet& subgrafo,
		VertexSet& clique, int cotaInferior, int cotaSuperior, int tolerancia)
{
	_stack[0]->_R->copyAll(subgrafo);
	_stack[0]->nodeidx = 0;
	_stack[0]->_R->nodeView(_stack[0]->nodeidx, *_stack[0]->Rnode);
	_stack[0]->_W->removeAll();

	int weightC = 0;
	int lastnode = NODEIDX(_n-1);

	int curC[_n];		//!< Current clique
	int topC;
	Doll * doll;
	int	v;
	int submax;

	//initialize with a maximal clique
	int expsz = 0;
	int i = expClique(*_stack[0]->_R, curC, &expsz, pesos);
	if (i > cotaInferior) {
		cotaInferior = i;
		clique.removeAll();
		for (int h = 0; h < expsz; h++)
			clique.add(curC[h]);

		int ub = greedyColoring(*_stack[0]->_R, pesos);
		if (cotaSuperior > cotaInferior + ub)
			cotaSuperior = cotaInferior + ub;

		if (cotaSuperior-cotaInferior <= tolerancia)
			goto TERMINATE;
	}

	topC = 0;
	doll = _stack[topC];
	v = doll->nextUpTo(lastnode);
	submax = 0;
	while (v >= 0 && submax + pesos[v] <= cotaInferior)
	{
		submax += pesos[v];
		_C[v] = submax;

		v = doll->nextUpTo(lastnode);
	}

	if (v < 0)
		goto TERMINATE;

	weightC = pesos[v];
	curC[topC++] = v;
	_stack[topC]->startUp(v, *_stack[topC-1]->_W);
	doll = _stack[topC];

	while (topC >= 0 && cotaSuperior-cotaInferior > tolerancia) {
		if (topC == 0)
			_C[curC[topC]] = submax;
		if (topC > 0 && doll->_R->isEmpty())
		{
			if (weightC > submax)
			{
				submax = weightC;

				if (submax > cotaInferior) {
					int expsz = topC;
					cotaInferior = expClique(*_stack[0]->_R, curC, &expsz, pesos);
					clique.removeAll();
					for (int h = 0; h < expsz; h++)
						clique.add(curC[h]);

					_stack[0]->_R->add(curC[0]);
					int ub = greedyColoring(*_stack[0]->_R, pesos);
					if (cotaSuperior > cotaInferior + ub)
						cotaSuperior = cotaInferior + ub;
					_stack[0]->_R->remove(curC[0]);
				}
			}

			topC--;
			if (topC >= 0) {
				weightC -= pesos[curC[topC]];
				doll = _stack[topC];
			}
		}
		else
		{
			v = doll->nextUpTo(lastnode);
			if (v < 0) {
				topC--;
				if (topC >= 0) {
					weightC -= pesos[curC[topC]];
					doll = _stack[topC];
				}
			}
			else if (topC == 0 || weightC + _C[v] > cotaInferior)
			{
				curC[topC++] = v;
				weightC += pesos[v];

				//Calculate the new candidate vector
				_stack[topC]->startUp(v, *doll->_W);
				doll = _stack[topC];
			}
			else if (weightC + _C[v] > submax)
				submax = weightC + _C[v];
		}
	}

	TERMINATE:

	return cotaInferior;
}

double OsterCliquer::greedyColoring(VertexSet& R, const DoubleValueArray& pesos) {
	double weight = 0.0, cw;
	int vidx;
	int nodeidx = 0;
	bool nemptcol = true;
	int lastnode = NODEIDX(_n-1);

	_WW.removeAll();

	while (nemptcol) {
		_RR.copyAll(R);
		nodeidx = 0;
		nemptcol = false;
		cw = 0.0;
		while (nodeidx <= lastnode) {
			_RR.nodeView(nodeidx, _nodeRR);
			R.nodeView(nodeidx, _nodeR);
			_WW.nodeView(nodeidx, _nodeW);
			vidx = _nodeRR.least();
			while (vidx >= 0) {
				int v = NELEM(nodeidx)+vidx;
				nemptcol = true;
				if (pesos[v] > cw)
					cw = pesos[v];
				// note that vidx is removed from _RR only if clqOrStab=false
				_RR.removeAll(_grafo.vecinos(v));
				_nodeRR.remove(vidx);
				_nodeR.remove(vidx);
				_nodeW.add(vidx);
				vidx = _nodeRR.least();
			}
			nodeidx++;
		}

		if (nemptcol)
			weight+=cw;
	}

	R.addAll(_WW);

	return weight;
}

double OsterCliquer::expClique(VertexSet& R, int * C, int * topC, const DoubleValueArray& pesos) {
	_RR.copyAll(R);

	double weight = 0.0;

	int t = *topC;
	while (--t >= 0) {
		// note that vidx only remains in _RR if clqOrStab=false
		_RR.retainAll(_grafo.vecinos(C[t]));
		weight+=pesos[C[t]];
	}

	if (_RR.isEmpty())
		return weight;

	int lastnode = NODEIDX(_n-1);
	int nodeidx = 0;
	int v;
	while (nodeidx <= lastnode) {
		_RR.nodeView(nodeidx, _nodeRR);
		int vidx = _nodeRR.least();
		while (vidx >= 0) {
			v = NELEM(nodeidx)+vidx;
			C[(*topC)++] = v;
			weight += pesos[v];
			// note that vidx only remains in _RR if clqOrStab=false
			_RR.retainAll(_grafo.vecinos(v));
			_nodeRR.remove(vidx);
			vidx = _nodeRR.least();
		}
		nodeidx++;
	}

	return weight;
}

double OsterCliquer::optimize(const DoubleValueArray& pesos, const VertexSet& subgrafo,
		VertexSet& clique,
		double cotaInferior, double cotaSuperior, double tolerancia)
{
	_stack[0]->_R->copyAll(subgrafo);
	_stack[0]->nodeidx = 0;
	_stack[0]->_R->nodeView(_stack[0]->nodeidx, *_stack[0]->Rnode);
	_stack[0]->_W->removeAll();

	double weightC = 0.0;
	int lastnode = NODEIDX(_n-1);

	int curC[_n];		//!< Current clique
	int topC;
	Doll * doll;
	int	v;
	double submax;

	//initialize with a maximal clique
	int expsz = 0;
	double i = expClique(*_stack[0]->_R, curC, &expsz, pesos);
	if (i > cotaInferior) {
		cotaInferior = i;
		clique.removeAll();
		for (int h = 0; h < expsz; h++)
			clique.add(curC[h]);

		double ub = greedyColoring(*_stack[0]->_R, pesos);
		if (cotaSuperior > cotaInferior + ub)
			cotaSuperior = cotaInferior + ub;

		if (cotaSuperior-cotaInferior <= tolerancia)
			goto TERMINATE;
	}

	topC = 0;
	doll = _stack[topC];
	v = doll->nextUpTo(lastnode);
	submax = 0.0;
	while (v >= 0 && submax + pesos[v] <= cotaInferior)
	{
		submax += pesos[v];
		_C[v] = submax;

		v = doll->nextUpTo(lastnode);
	}

	if (v < 0)
		goto TERMINATE;

	weightC = pesos[v];
	curC[topC++] = v;
	_stack[topC]->startUp(v, *_stack[topC-1]->_W);
	doll = _stack[topC];

	while (topC >= 0 && cotaSuperior-cotaInferior > tolerancia) {
		if (topC == 0)
			_C[curC[topC]] = submax;
		if (topC > 0 && doll->_R->isEmpty())
		{
			if (weightC > submax)
			{
				submax = weightC;

				if (submax > cotaInferior) {
					int expsz = topC;
					cotaInferior = expClique(*_stack[0]->_R, curC, &expsz, pesos);
					clique.removeAll();
					for (int h = 0; h < expsz; h++)
						clique.add(curC[h]);

					_stack[0]->_R->add(curC[0]);
					double ub = greedyColoring(*_stack[0]->_R, pesos);
					if (cotaSuperior > cotaInferior + ub)
						cotaSuperior = cotaInferior + ub;
					_stack[0]->_R->remove(curC[0]);
				}
			}

			topC--;
			if (topC >= 0) {
				weightC -= pesos[curC[topC]];
				doll = _stack[topC];
			}
		}
		else
		{
			v = doll->nextUpTo(lastnode);
			if (v < 0) {
				topC--;
				if (topC >= 0) {
					weightC -= pesos[curC[topC]];
					doll = _stack[topC];
				}
			}
			else if (topC == 0 || weightC + _C[v] > cotaInferior)
			{
				curC[topC++] = v;
				weightC += pesos[v];

				//Calculate the new candidate vector
				_stack[topC]->startUp(v, *doll->_W);
				doll = _stack[topC];
			}
			else if (weightC + _C[v] > submax)
				submax = weightC + _C[v];
		}
	}

	TERMINATE:

	return cotaInferior;
}

